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Abstract. In the two-dimensional Ising model weak random surface field is predicted 
to be a marginally irrelevant perturbation at the critical point. We study this question 
by extensive Monte Carlo simulations for various strength of disorder. The calculated 
effective (temperature or size dependent) critical exponents fit with the field-theoretical 
results and can be interpreted in terms of the predicted logarithmic corrections to the 
pure system's critical behaviour. 

PACS numbers: 05.70Jk, 05.50.+q 
1. Introduction 

In an inhomogeneous system the local critical behaviour near localized or extended 
defects may differ considerably from the bulk critical behaviour in the regular lattice 
(for a review, see One possible source of inhomogeneity is quenched (i.e., 

time-independent) randomness, which can be localized at the surface of the system 
(fluctuating surface coupling constants [2-6], microscopic terraces at the surface [7]) or 
at a grain boundary in the bulk of the system. It is known experimentally [8-10] that 
impurities may diffuse from inside the sample and segregate on the surface or at grain 
boundaries. In adsorbed systems, quenched disorder is naturally present since adatoms 
may bind randomly on equivalent surface sites 

In a theoretical description of the local critical behaviour of these systems, close 
to the bulk critical point, one can use a coarse-grained picture in which quenched 
randomness couples to some local operator. The local operator considered in this paper 
is the surface order parameter, hence the perturbation is described by the introduction 
of random fields (RFs) localized at the surface. Usually the RF has zero mean and its 
variance is used to characterize the strength of disorder. 
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In the weak-disorder limit, the relevance or irrelevance of the perturbation can be 
analyzed by making use of a Harris-type criterion [T^j . The condition for the irrelevance 
of RFs on a defect with dimension d — 1 can be expressed in terms of the decay exponent 
for the local order parameter correlations in the pure system as 



A plane of RFs in the bulk often constitutes a relevant perturbation as it is the 
case for the two-dimensional (2D) Ising model with rj = 1/4. Thus a new fixed point 
appears, which controls the local critical behaviour. This fixed point is expected to be 
a surface one since RFs tend to destroy the local order and the bulk defect then acts as 
an effective cut. 

Surface RFs are irrelevant for the 3D Ising model as noted and demonstrated 
through Monte Carlo (MC) simulations ^H]- Curiously, in the case of a system with 
continuous symmetry, like the 3D Heisenberg model, surface RFs destroy the bulk long- 
range order jT3] according to Imry-Ma arguments ^21, although the perturbation is 
irrelevant at the ordinary surface transition according to (^. Among 2D systems the 
Ising model represents the borderline case, since ri\\ = 1 [T^. For this model, field- 
theoretical investigations fUEI predict that weak surface RF is a marginally irrelevant 
perturbation. Consequently, the surface critical properties of the random model are 
characterized by the critical singularities of the pure model supplemented by logarithmic 
corrections to scaling. 

These theoretical predictions have not yet been confronted with the results of 
numerical calculations. In general, the observation and characterization of logarithmic 
corrections to scaling by numerical methods are notoriously difficult tasks, particularly 
in systems with quenched disorder. In this respect, a well-known example is the diluted 
2D ferromagnetic Ising model, for which the accurate form of the singularities was long 
debated [18-23]. 

In this paper we present the results of a numerical study of the surface critical 
behaviour of the 2D Ising model in the presence of random surface fields. In 
section 2 we present the model and the known results about its critical properties. 
Then, through intensive Monte Carlo simulations, we determine effective surface 
magnetization exponents in two different ways. In section 3, they are obtained as a 
function of the deviation from the critical temperature. In section 4, we use a small 
homogeneous surface field at the critical point to deduce size- dependent exponents from 
the magnetization profiles. In section 5 we discuss the agreement between theoretical 
and numerical results. 

2. The model and its predicted surface critical properties 

We consider the Ising model on a L x M square lattice with the Hamiltonian 
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^ , \ hs + h with probability p=l/2 
*' ' * ] hs — h with probability p=l/2 

where Sij = ±1. J is the first-neighbour exchange interaction. The RF hs ± h acts 
on the surface spins in the columns at i = 1 and i = L and periodic boundary 
conditions are used in the vertical direction. Our main interest is to calculate the 
averaged magnetization per column, mj = (| J2j 

For the pure system, i.e. with vanishing surface field, the surface magnetization, 
m-s = mi = mi is exactly known in the thermodynamic limit (L, M —>■ oo) jl^: 

sinh(2ir) - 



^s,purc 



coth(2ir) 



(3) 



cosh(2/s:) - 1 

in terms of i^' = J/k^T, where T is the temperature. At the critical point with 
sinh(2frc) = 1 the surface magnetization vanishes as 

"^s.pure ~ mot^^'^ , (4) 

in terms of the reduced temperature, t = (Tc — T)/Tc, and ml = 4(v^ + l)Kcla.{2Kc). 
The relevant length scale is the bulk correlation length which, for t > 0, is given by [21] 

e=[21n(sinh(2ir))]~^ , (5) 

with the lattice constant for unit length. The correlation length diverges at the critical 
point as ^ ~ [2^/2Kct]~^. Thus the reduced temperature and the length scale are related 
by 

-Int ~ .913 + ln^ . (6) 

In the presence of random surface fields there are no exact results available. In this case 
one can use the replica trick to transform the semi-infinite system with a random surface 
field into n semi-infinite replicas, coupled two-by-two through their surface spins, via 
nearest-neighbour interactions proportional to h"^. The average properties of the random 
system are obtained in the limit n ^ 0. 

The surface critical properties of the system have been studied via two different 
methods, both using the differential renormalization group (RG) techniques where the 
lengths are rescaled by a factor e' . The surface coupling between the replicas transforms 
as 

In the first approach ^7] a conformal mapping is used at the bulk critical point, 
with /is = 0, to transform the n semi-infinite replicas into n infinite strips with width 
L, which are coupled to each other at both surfaces through h"^. The behaviour of 
the inverse correlation length is studied using degenerate perturbation theory to second 
order in h"^. From the transformation of the inverse correlation length on the strips 
under rescaling by e' = L, with n = 0, and using the gap-exponent relation jSHj, one 
can identify the L-dependent, effective decay exponent t]\\ which is given by 

'(,l = l + i^. (8) 
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to leading logarithmic order. 

In the second approach ^1] the behaviour under rescaling of the homogeneous part 
of the surface field is determined as 

The surface free energy density transforms as 

Ut, = hs{l), h^l)] , (10) 

where u is the correlation length exponent. Using (j^I), the surface magnetization reads 



m,{t,h ) = — 
oh 



p-'/2 

With e' = ^ ~ t^'^ and u = 1 according to (0), ignoring higher order corrections, one 
obtains 

Thus the critical singularity of the pure model is supplemented by a logarithmic 
corrections to scaling. From a practical point of view, one can define temperature- 
dependent effective exponents through 

lnK(t(l + ^))/ms(t(l-(5))] _ 1 / 1 \ 
ln[(l + 5)/(l-5)] -2i^^ + rM+"V ' 

with 5 — > 0. The last expression gives the leading logarithmic correction following from 
equation (fT^ . 

Taking into account the scaling relation r^y = 2(3^/^ with v = 1, the effective 
exponents in equations (fT^ and (jS)) correspond in terms of the relevant length scales, 
L~^i^ l/t. 

The following two sections are devoted to a numerical test of the validity of these 
theoretical results. 



3. Effective exponents 

The surface critical exponents are deduced from the temperature dependence of the 
magnetization rui in the surface layers (for a review, see [26j). We set the homogeneous 
surface field to zero, /ig = 0, the strength of the random surface field ranging from 
h = 0.6 to h = 1.5, and take a finite reduced temperature, t > 0. Systems of square and 
rectangular shapes containing L x M spins, with L and M ranging from 50 to 1000, 
have been studied using the standard single-spin-flip method. Although systems with 
rectangular shapes (M < L) lead to reduced finite-size effects, the fraction of surface 
spins is smaller than for a square system and more runs are needed to achieve the same 
accuracy for the surface magnetization which is self-averaging. Thus we worked with 
square systems to spare computer time. The final data are obtained after averaging 
over at least 1000 different runs with different realizations of the random surface field. 



Two-dimensional Ising model in a random surface field 



5 



For every run, time average has been taken over a few 10^ Monte Carlo steps per spin 
after equilibration. 

As an illustration we present rui aX t = 0.05 and t = 0.02 in figure Q for different 
strengths of the RF. The profile of the pure system is shown for comparison. For a 
given t and different values of h, mi displays a plateau around i = L/2 for large enough 
systems. It corresponds to the bulk magnetization, mt, since its height is independent 
of h. If we approach the surfaces close enough, i,L — i < ^, we enter in the surface 
region where the value of the magnetization is rapidly decreasing to its surface value, 
mg. As seen in figure ^ for a given t the surface magnetization mg and the inverse size 
of the surface region are decreasing with increasing disorder strength h. 




Figure 1. Magnetization profiles with random surface fields h — 0, 0.5 and 1 at 
t = 0.05 and t — 0.02. The data have been obtained for a system with 300 x 300 spins. 

According to finite-size scaling, in a large but finite system, sufficiently close to its 
critical point, ms(t) behaves as (— t)^"/(^/L), where the scaling function, /(x), tends 
to a constant for small values of its argument. In the actual calculations, we approach 
the transition point only to such a distance that the finite-size effects remain negligible 
and effective surface exponents are calculated using (fT!^ . In practice, finite-size effects 
have been circumvented by adjusting the size of the sample in a standard approach 
For a given value of t, data obtained for different system sizes are compared. Away 
from the critical point these data agree as long as the correlation length is less than 
the extent of the smaller system. Closer to Tc the correlation length increases and at 
some stage it gets comparable to the size of the smaller system. Finite-size effects then 
show up by a characteristic fast drop of the effective exponent [Ij. The smaller system 
is then discarded and the procedure is continued with two system sizes which still yield 
identical data at that temperature. This approach is somehow cumbersome but assures 
that the final data are essentially free of finite-size effects. 
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The effective exponents are shown in figure |21 for different values of the strength of 
the RF. Here, in order to check the form of the logarithmic corrections in equations (jHl) 
and (|T^ . Ps is plotted as a function of 1/| lnt|. For a given t the effective exponents are 
increasing with h as expected since RFs decrease local order. On the other hand, for 
a given h, Psit) first shows a monotonic increase when t decreases and its value passes 
over the pure system's surface exponent, jSg = 1/2. Then, by further decreasing the 
reduced temperature, /3s (t) seems to approach a maximum value. This saturation effect 
is more evident for large values of h. Unfortunately, the size limitation did not allow 
us to approach the transition point close enough to follow the predicted decrease of the 
effective exponent. 

Here, in order to compare the numerical results with the theoretical predictions 
and to extrapolate our data to t — 0, we use the following expression, 

= + V' • 

which contains the leading analytic correction which follows from equation Q and the 
leading logarithmic corrections to the fixed-point singularity given by (fT^ . Note that 
the two corrections have different signs and their competing effect results in the non- 
monotonic temperature dependence of the effective exponent, /3s(t). For a given RF, we 
have fitted the surface magnetization data to the form given in (fT^ with the amplitudes 
a and b as free parameters. From this, the effective surface magnetization exponent was 
calculated and used to extrapolate the data points in figure El The data extrapolate to 
the value of the pure system jSg^pure = 1/2 after an "overshooting effect". A quite similar 
tendency was observed in [20] for the bulk magnetization exponent of the 2D random 
bond Ising model. In our calculations, however, the maximum value is almost reached 
for h > 1, which was not possible for the random bond model. 

In order to check the leading logarithmic correction to the effective exponent, we 
have calculated the difference of /3s (t) and its value in the pure system, /3s,pure(^), as 
calculated from via equation (|T3|l . This difference, which is plotted in the inset of 
figure 121 no longer contains the leading analytic correction; therefore we expect to be 
able to compare it with the theoretical prediction in equation ()13|). As seen in the inset 
for the random surface field, h = 1, the corrections are compatible with theory, although 
much larger systems are needed in order to reach the asymptotic regime. 

The form of the logarithmic corrections has been analyzed in still another way 
by forming the ratio, r{t,h) = ms,pure/"^s5 of the surface magnetizations in the pure 
and in the disordered systems. In this way the leading analytic correction to scaling 
is eliminated. As shown in figure El the square of r(t, h) has an asymptotic linear 
dependence on Int, the slope of which is proportional to /i^, as shown in the inset of 
figure El This is in complete agreement with the theoretical prediction in equation (fT^ . 
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1/|lnt| 

Figure 2. Effective surface exponent for different strengths of the random surface 
field. The two first points (grey circles) for h — 1 were deduced from the short-distance 
behaviour of the critical profiles discussed in section 4. The broken lines fitting the 
numerical data correspond to the formula in equation H14() . with a — —0.60, b ~ —0.23 
{h = 0.6), a = -0.53, b = -0.40 {h = 0.8), a = -0.40, b = -0.47 {h = 1.0) and 
a = —0.47, b ~ —0.72 {h = 1.5). The inset gives the difference between the effective 
exponents for the random {h — 1) and the pure systems. The straight line gives the 
predicted leading logarithmic correction. 
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Figure 3. Square of the ratio r{t, h) of the surface magnetizations in the pure and in 
the disordered Ising model as a function of Int. The slopes are proportional to ft.^ as 
shown in the inset. Error bars are much smaller than the symbol sizes. 
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4. Critical profiles 

In this section we study the system at the critical point, t = 0; however, in the presence 
of a small homogeneous surface field, hs <^ h. A typical magnetization profile in the 
system with M <^ L is shown in the inset of figure IH These data have been obtained 
with the Swendsen-Wang algorithm with a layer of ghost spins next to the surface |27j . 
At least 1000 runs with different realizations of the random surface field were performed, 
the time average of every run resulting from typically 3 ■ 10^ MC updates. 

As known from an analysis of the non-random system j2Zj the profile first increases 
close to the surface and then decreases in the bulk. The surface critical exponent, 
f3s, influences the form of the initial part and can be extracted from it. The non- 
vanishing surface field introduces a new surface length scale. Is, which in 2D is given 
by Is ~ hs ^^^^ and scales as ~ for the Ising model. 
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Figure 4. Ratio of the initial part of the critical magnetization profiles calculated in 
the random and in the pure system in the presence of a small homogeneous magnetic 
field. The initial slope in the log-log plot corresponds to the difference in the effective 
magnetization exponents (see the text). Inset: one half of the magnetization profile in 
the random system. 

The initial part of the profile is restricted to i -C M, /g and, according to finite-size 
scaling theory it behaves as 

m,{M,ls)=t-^-g{t/M,t/ls) , (15) 

where Xin = /?/z/=l/8is the bulk magnetization scaling dimension. For small values of 
its arguments the scaling function g is expected to factorize as goii/ M)gs{i/ls)- For the 
pure 2D Ising model the second term is logarithmic j^, gsii/h) ~ ~ lii(^^s)' 

which makes the numerical analysis difficult. For the Ising model with RF on the surface, 
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the first term is expected to behave as 

This form incorporates logarithmic corrections which, with i = 1 and M ~ ^, are in 
agreement with the form of the surface magnetization given in equation (|T^. Now, in 
analogy with (|T3|l . an effective surface magnetization exponent can be defined as 

_ H9o{ya-5))/goiyii + S)) 

V ln[(l + (5)/(l-5)] ^ ^ 

which, finally when ?/ <^ 1, will be a function /3s (M) of the characteristic length of the 
problem. 

In order to get rid of the logarithmic factor, gj^ijl^^, we have calculated the ratio 
mi(M)/mj_pure(M) of the initial profiles in the random and pure systems. Assuming 
that g^ has the same logarithmic singularity for both, we arrive to the conclusion that 
the ratio of profiles scales with the difference of the effective exponents, A/3s(M) = 

/?s(M)-/?s,pure(M). 

In the actual calculation we set = .01 and h = 1, and performed MC simulations 
on systems with sizes M = 16 and 64 and different aspect ratios a = L/M. For 
a = 4 and 8, the initial part of the profiles turned out to be indistinguishable. The 
ratios calculated for the largest a are presented in figure El and, from the extrapolated 
initial slopes in a log- log plot, the differences of the effective exponents take the values 
A/?s(16) = 0.08(1) and A/3s(64) = 0.07(1). In order to compare these estimates with the 
effective exponents obtained in section 3 with a finite t, we use the correspondence of 
equation with M = ^. The data points obtained in this way are inserted in figure El 
They seem to fit very well with the predicted theoretical curve and are located in the 
descending part of the curve. Therefore this calculation gives further support to the 
theoretical results about the form of the logarithmic corrections. 



5. Discussion 



Marginally irrelevant operators are responsible for logarithmic corrections to scaling, the 
form of which can be often predicted by field theory and conformal invariance. According 
to second-order perturbation theory, the random surface field is expected to be such a 
marginally irrelevant operator at the surface fixed point of the 2D Ising model. This 
conjecture, which is made on the basis of the replica trick and in the weak disorder limit, 
is confronted here with the results of extensive MC simulations for varying strength of 
the disorder. The calculated effective surface magnetization exponent, which depends 
either on the distance t from the critical point or on the finite size M of the critical 
system, varies with these parameters. Since this variation is non-monotonic, a direct 
extrapolation to the fixed-point values cannot be made from the data available on finite 
systems. However, the variation of the effective exponents is in good agreement with the 
theoretical form, which contains the predicted logarithmic correction to scaling to the 
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pure systems critical behaviour. For the largest systems and for the strongest random 
fields, the numerical results are not too far from that obtained perturbatively in linear 
order. Therefore we interpret our results as numerical evidence in favour of the validity 
of the field-theoretical predictions. 
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